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Motivated by recent experiments in ultracold atomic gases that explore the nonequilibrium dy- 
namics of interacting quantum many-body systems, we investigate the opposite limit of Landau's 
Fermi liquid paradigm: We study a Hubbard model with a sudden interaction quench, that is the in- 
teraction is switched on at time t = 0. Using the flow equation method, we are able to study the real 
time dynamics for weak interaction U in a systematic expansion and find three clearly separated time 
regimes: i) An initial buildup of correlations where the quasiparticles are formed, ii) An intermediate 
quasi-steady regime resembling a zero temperature Fermi liquid with a nonequilibrium quasiparticle 
distribution function, iii) The long time limit described by a quantum Boltzmann equation leading 
to thermalization of the momentum distribution function with a temperature T oc U. 



The investigation of interacting quantum many- 
particle systems in nonequilibrium has recently attracted 
a lot of attention. A simple way to excite a system from 
its ground state is an interaction quench, a sudden switch 
of parameters in the Hamiltonian. The time evolution 
of the initial state is then generated by the quenched 
Hamiltonian, for which the initial state is generically not 
an eigenstate. Recent experiments have implemented 
quenches of ultracold atoms loaded on optical lattices 
and observed remarkable subsequent dynamics described 
as iterated 'collapse and revival' of the initial superfluid 
phase [TJ Ej. Yet their theoretical description remains 
a challenge since many well-established equilibrium the- 
oretical methods fail in nonequilibrium. From a theo- 
retical point of view, the long-time limit poses partic- 
ularly intriguing questions: Will an interacting closed 
quantum system prepared in some generic initial state 
equilibrate, that is behave like the equilibrium system 
with some nonzero temperature after waiting sufficiently 
long? In nonlinear classical systems similar questions 
have been addressed in a multitude of publications since 
the seminal work by Fermi, Pasta and Ulam [3]. Non- 
equilibration has been linked to integrability since an in- 
tegrable system is constrained by an infinite number of 
conservation laws. 

However, much less is known about quantum systems. 
Since a pure state remains a pure state under unitary 
time evolution, the concept of thermalization is only 
meaningful for suitable observables. First theoretical re- 
sults have shown that observables may approach limiting 
values or exhibit persistent oscillations which, even when 
time-averaged, do not match with equilibrium properties 
[HE]. A proposition by Rigol et al. [6] gave a statisti- 
cal description for the stationary state of an integrablc 
system in terms of a generalized Gibbs ensemble. Con- 
ditions for the applicability or non-applicability of this 
scenario have been clarified in [7] and specific results 
have been obtained for the Luttinger model [8 , hard core 
bosons in one dimension [9j [10] and the infinite dimen- 
sional Falicov-Kimball model [11 . While the concept of a 



generalized statistical ensemble proved helpful even for a 
less restrictive set of constraints [H] , the role of integra- 
bility has been questioned by further numerical works: 
Breaking the integrability of spinless fermions on a Id 
lattice has not altered relaxation to a non-thermal state 
[12] . Similarly, for the non-integrable Id Bose- Hubbard 
model signatures of thermalization could only be found 
for a limited regime of quenches, while others seemed to 
drive the system to non-thermal stationary states |13j . 
Exact results have been obtained for the opposite case 
of quenches from the Mott phase to the noninteracting 
Hamiltonian and show relaxation of local observables to 
a nonequilibrium steady state |14j . 

Motivated by these questions, we study an interaction 
quench in a Fermi liquid in d > 1 spatial dimensions, that 
is we suddenly switch on the interaction at time t = 0. 
This is the extreme opposite limit of Landau's adiabatic 
switching on procedure, where one finds the celebrated 
one to one mapping between physical electrons and quasi- 
particles. In the sudden quench scenario, the system is 
prepared as the zero temperature ground state of the 
noninteracting Fermi gas at times t < 0, and then, for 
t > 0, subject to the time evolution with respect to the 
interacting Hamiltonian. We find three regimes of the 
time evolution, that are well separated for weak interac- 
tion: An initial quasiparticle formation regime, followed 
by a quasi-steady intermediate regime resembling a zero 
temperature Fermi liquid, and a long-time thermalization 
regime where the momentum distribution function equili- 
brates. Concretely, we investigate the fermionic Hubbard 
model at half filling described by the following Hamilto- 
nian (Fermi energy ep = 0) 

H(t)= ]T e k :cl^ ka :+Q(t)UY,^-\){n ii - 1 -) 

(1) 

and work out the time-dependent momentum distribu- 
tion functions Nk(t). Notice that this system is clearly 
non-integrable for d > 1 and one therefore expects 
generic behavior. Most of our results are obtained in 
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FIG. 1: The Heisenberg equation of motion for an observable 
O is solved by transforming to the B = oo eigenbasis of the in- 
teracting Hamiltonian H (forward transformation), where the 
time evolution can be computed easily. Time evolution intro- 
duces phase shifts, and therefore the form of the observable 
in the initial basis B — (after a backward transformation) 
changes as a function of time. 



the limit of high dimensions [15 , but the calculation also 
applies to finite dimensions with the same conclusions up 
to quantitative details. 

We study the above real time evolution problem by 
using the approach introduced in [16]. One solves the 
Heisenberg equations of motion for the operators that one 
is interested in by performing a unitary transformation 
to an (approximate) eigenbasis of the interacting Hamil- 
tonian. There one can easily work out the time evolution 
and then transform back to the original basis where the 
initial state is specified. In this manner one induces a 
solution of the Heisenberg equations of motion for an op- 
erator in the original basis but without secular terms, 
which are usually a major problem in other approxima- 
tion schemes [T7]. Fig. [I] gives a sketch of our approach. 
Notice that the same general idea was recently also used 
by Cazalilla to study the behavior of the exactly solvable 
one-dimensional Luttinger model subject to a quench [5]. 

Since our model is non-integrable, we implement the 
above diagonalizing transformation by the flow equa- 
tion method [HJ [19], which permits a systematic con- 
trolled expansion for many equilibrium and noncqui- 
librium quantum many-body problems |19j . One uses 
a continuous sequence of infinitesimal unitary transfor- 
mations parametrized by a parameter B with dimen- 
sion (energy) -2 that connects the eigenbasis of the free 
Hamiltonian (B = 0) with the energy diagonal basis of 
the interacting Hamiltonian (B — oo). Each infinites- 
imal step of the unitary transformation is defined by 
the canonical generator rj(B) = [Hq(B), Hi nt (B)], where 
Hq(B) is the diagonal and H lnt {B) the interacting part 
of the Hamiltonian. This generator rj(B) has the re- 
quired property of making H(B) increasingly energy di- 
agonal for B — > oo [18J. All operators O(B) (including 
the Hamiltonian itself) flow according to the differen- 
tial equation dO(B)/dB = [rj(B),0(B)]. Higher order 
terms generated by the commutator are truncated af- 
ter normal-ordering (denoted by : :) and the flow equa- 
tions decompose into a set of ordinary differential equa- 
tions resembling scaling equations in a renormalization 



approach. However, contrary to conventional renormal- 
ization schemes which reduce the size of the effective 
Hilbert space, the flow equation approach retains the full 
Hilbert space, which makes it particularly appropriate 
for nonequilibrium problems (for more details sec [19 ). 

Flow equations for the Hubbard model. First we work 
out the diagonalizing flow equation transformation for 
the Hubbard Hamiltonian. The expansion parameter is 
the (small) interaction U and normal-ordering is with 
respect to the zero temperature Fermi-Dirac distribution: 



H(B) = ]T e k :clc ka : 

fc<r=T4 

+ E U p'p<i'i( B ) 4, 

p'pq'q 



(2) 



T c pT c «'1 c ?1 : 



with U p i pq 'q(B = 0) = U. The flow of the one-particle en- 
ergies and the generation of higher normal-ordered terms 
in the Hamiltonian can be neglected since we are in- 
terested in results in second order in U. The flow of 
the interaction is to leading order given by U p i pq ' q {B) — 

C/exp(— BA plpq , ) with an energy difference A p > pq i q = f 

£p' — £p + Eg' — £q- 

Next we work out the flow equation transformation for 
the number operator N k -\(B) = C k AB)C k AB), which 
can be obtained from the transformation of a single cre- 
ation operator C\^(B). Under the sequence of unitary 
transformations the operator changes its form to describe 
dressing by electron-hole pairs. A truncated ansatz reads: 



Cj T (B) = ft fc (B) C t 



p'i'v 



p'+q' :C p'1 C q'l C pi : 



. (3) 

We introduce the zero temperature momentum distri- 
bution function of a free Fermi gas n k , define n k d = 

n P ' n q' U p + 



def 



1 — n k and a phase space factor Q p > pq i [n] 
n p ,n q ,n p . The flow equations for the creation operator 



d -^> =UJ2 M k p , q , p {B) A kp , pq , e - BA W Qplp , r[lll 



p'q'p 



dM k plq , p {B) 
dB 



h k (B)UA p , pq , k e- BA *'»«"° (4) 



Here and in the ansatz (|3J) we have only taken into ac- 
count the terms that are required to describe the momen- 
tum distribution function up to second order in U. The 
initial conditions for the above transformation of are 

MO) = 1 and M* qlp (0) = (i.e., C^B = 0) = 4 T ), 
and we denote the asymptotic values from the solution 
of Q by h k (B = 00, t = 0) and Mp q , p (B = 00, < = 0). 
Time evolution according to Fig. [TJyields h k (B = 00, t) — 
h k {B = 00, t = 0)e- lCkt and M^, q)p (B = 00, t) = 
M P 'q' P ( B = = 0)e -i( V+V- £ p)*, which are then 

input as the initial conditions of the system of equations 
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Q at B = oo. Integrating back to _B = gives the 
time evolved creation operator in the original basis, and 
it is straightforward to evaluate the time dependent mo- 
mentum distribution function with respect to the initial 
Fermi gas state [20] . 

Nonequilibrium momentum distribution function. One 
finds the following time-dependent additional term to the 
distribution n k of the free Fermi gas in 0(U 2 ): 

AN^(t) = Nr Q (t)~n k (5) 
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The phase space factor J k (E; n) resembles the quasipar- 
ticle collision integral of a quantum Boltzmann equation: 
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For computational convenience we use the limit of infi- 
nite dimensions, specifically a Gaussian density of states 
p(e) = exp(-(e/t*) 2 /2) /V2nt* [15]. In the sequel 
pp = p(e — 0) denotes the density of states at the Fermi 
level. Results from a numerical evaluation of the above 
scheme for three time steps are presented in Fig. [2] 

Equilibrium momentum distribution function. Eqs. Q 
can also be used to evaluate the equilibrium distribution 
function, which will later be important for comparison. 
In fact, the asymptotic value h kp (B = oo) at the Fermi 
energy is directly related to the quasiparticle residue (Z- 
factor), Z EQU = [h kF (B = oo)] 2 [19 . It is easy to solve 
Q analytically at the Fermi energy for zero temperature 
in 0(U 2 ) and one finds for momenta k infinitesimally 
above or below the Fermi surface 



AN, 



EQU _ 



= -U d 



dE 



■ME] n) 
{tk-E) 2 



(6) 



consistent with a conventional perturbativc evaluation. 

Short-time correlation buildup. The numerical evalua- 
tion of the momentum distribution function depicted in 
Fig.[2]shows the initial buildup of a correlated state from 
the Fermi gas. For times < t < p~j? l U~ 2 one observes 
a fast reduction of the Fermi surface discontinuity and 
1/t oscillations in the momentum distribution function. 
This short time regime can be understood as the forma- 
tion of quasiparticles from the free electrons of the initial 
nonintcracting Fermi gas. 

Intermediate quasi-steady regime. For times t of or- 
der pp 1 !/^ 2 the sinusoidal time dependence in |5j) gen- 
erates an increasing localization in energy space, which 
eventually becomes a (5-function (Fermi's golden rule). 
There are no further changes in the momentum distribu- 
tion function for times t > Pp^U~ 2 in the present order 
of the calculation. For momenta k infinitesimally above 
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FIG. 2: (a)-(d): Time evolution of iV NEQ (e) plotted around 
the Fermi energy for pfU = 0.6. A fast reduction of the dis- 
continuity and 1/t-oscillations can be observed. The arrow in 
(d) indicates the size of the quasiparticle residue in the quasi- 
steady regime. In (e) the universal curves for AiVfc = Nk — nu 
are given for both equilibrium and for the nonequilibrium 
quasi-steady state in the weak-coupling limit. 



or below the Fermi surface one then finds from (pil): 



Ai\£ EQ (i-oo) - -W 2 



dE 



1_ J k {E; n) 
2 (e k -E) 2 



2 A Af QU 



(7) 



since sin 2 in (jE]) yields a factor 1 /2 in the long time limit. 
In the quasi-steady state the momentum distribution 
function is therefore that of a zero temperature Fermi 
liquid. However, from ^ one deduces that its Z-factor 
is smaller than in equilibrium, 

l _ Z NBQ = 2 (1-Z E Q U ). 

This factor 2 implies a quasiparticle distribution function 
in the vicinity of the Fermi surface in the quasi-steady 
state equal to the equilibrium distribution function of 
the physical electrons, _/V^ P NE( ^ = _/V ECJU , as opposed to 
its equilibrium distribution, 7\rQ p E Q u = 0(fc F — k). 

Remarkably, Cazalilla's findings [S] for the interaction 
quench in the Luttinger model mirror these features: the 
critical exponent describing the asymptotic behavior of 
the electronic Green's function differs from the equilib- 
rium result. As Cazalilla points out this corresponds to a 
non-equilibrium distribution for the bosonic modes after 
bosonization. A main difference between the Luttinger 
liquid and the Fermi liquid cases follows from the inte- 
grability of the Luttinger liquid with an infinite number 
of conservation laws, which make this regime stable for 
t — * oo. For the Fermi liquid, on the other hand, on-shell 
interactions lead to thermalization as we will see next. 

Thermalization. The previous flow equation calcula- 
tion of the real time dynamics contains all contributions 
to the time evolution for times smaller than pp 3 U~ 4 . For 
the long time dynamics one generally expects a quantum 
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Boltzmann equation (QBE) to be a valid description [2"T] 

= - PF U 2 J k (E = e k ,N^(t)) . (8) 

Here the quasiparticle momentum distribution function 
-W* P ' NE derived above serves as the initial condition. 
Because jv^ p ne( ^ allows nonzero phase space for scat- 
tering processes in the vicinity of the Fermi surface (orig- 
inating, ultimately, from the factor 2 in ([7])), the initial 
quasiparticle distribution function starts to evolve on the 
time scale t oc p F 3 U~ 4 . This implies that the quasi- 
steady electron distribution function depicted in Fig. [2ji 
starts to decay on this time scale and one approaches 
a Fermi-Dirac distribution (being the only stable fixed 
point of ([8])) with a nonzero temperature T. 

The above scenario fits well into the picture of noncqui- 
librium field theories describing, e.g., the early universe 
[IT] . The excitation energy of the initial quantum state 
(the Fermi gas) with respect to the equilibrium ground 
state of (JlJ is E cyi = appU 2 in the weak interaction limit 
with some lattice-dependent constant a > 0. The short- 
time correlation buildup corresponds to prethermaliza- 
tion, where kinetic and interaction energy in flow fr om 
to E™f Q = -2a PF U 2 and E^ Q = 2ap F U 2 . This fol- 
lows immediately from the Feynman-Hcllman theorem 
and the fact that the total energy remains zero for all 
times. -Ejnf equals the equilibrium interaction energy, 
while -E£! EQ = + E cx . Kinetic and interaction en- 

ergy then remain constant throughout the quasi-steady 
regime and the long-time limit, and therefore the system 
has prethermalized for these average quantities. In the 
thermalization regime the system redistributes its addi- 
tional excitation energy E cx in the kinetic energy over 
the different momenta and reaches a Fermi-Dirac distri- 
bution with temperature T cx U. 

Higher order flow equations. Clearly, it would be desir- 
able to derive ^ within the framework of the real time 
flow equation calculation. However, a calculation to or- 
der U 4 is beyond the scope of the present work. Still, 
one can identify a particular contribution in fourth order 
leading to a finite lifetime of order p F 3 U~ 4 for an elec- 
tron at the Fermi surface, which is consistent with the 
dynamics implied by the QBE. The short time evolution 
of the system for times smaller than p F 3 U~ 4 obtained 
from the full solution of the Heisenberg equations of mo- 
tion therefore matches the long time dynamics described 
by the QBE, and we have a consistent picture on all time 
scales. Another effect of the fourth order contributions 
is that the sharp Fermi edge of the quasi-steady state 
gets smeared out on an energy scale ppU 4 , which, how- 
ever, does not essentially modify our previous conclu- 
sions. Therefore, strictly speaking the discontinuity of 
the momentum distribution function disappears immedi- 
ately for t > 0, but this effect only becomes noticeable 
for times of order p F 3 U~ 4 . 



Conclusions. We have discussed the real time evo- 
lution of the Hubbard model with a sudden interac- 
tion quench for a weak interaction U. Ultimately, the 
system completely thermalizes its excitation energy E ex 
and reaches a temperature T oc U. This thermaliza- 
tion regime only sets in on the time scale p F 3 U~ 4 . This 
follows from the observation that the short time behav- 
ior up to times of order p F l U~ 2 amounts to quasipar- 
ticle formation with a momentum distribution function 
with a discontinuity at the Fermi energy. Therefore, 
a quasi-steady prethermalized state emerges for times 
Pp 1 !/^ 2 < t < p F 3 U~ 4 . Its momentum distribution 
function looks like a zero temperature Fermi liquid, but 
with the wrong quasiparticle residue with respect to the 
interacting ground state. It is this noncquilibrium quasi- 
particle residue that allows for phase space for scattering 
processes in a quantum Boltzmann equation description 
for times t > p F 3 U~ 4 , which then leads to thermalization 
of the momentum distribution function. 
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